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Abstract. 

Spurred by the improved measurement sensitivity in the very-high-energy (VHE: >100 GeV) 7-ray band, we assess the fea- 
sibility of detection of the nearby starburst galaxy M 82. VHE emission is expected to be predominantly from the decay of 
neutral pions which are produced in energetic proton interactions with ambient protons. An estimate of VHE emission from 
this process is obtained by an approximate, semi-quantitative calculation, and also by a detailed numerical treatment based on a 
convection-diffusion model for energetic electron and proton propagation and energy losses. All relevant hadronic and leptonic 
processes are considered, gauged by the measured synchrotron radio emission from the inner disk region. We estimate an inte- 
grated flux /(>100GeV)~2x 10 -12 cm _2 s~ 1 , possibly detectable by the current northern-hemisphere imaging air Cherenkov 
telescopes, MAGIC and VERITAS, and a good candidate for detection with the upcoming MAGIC II telescope. We also esti- 
mate /(> 100MeV)~10 -8 cm _2 s _1 , a level of emission that can be detected by GLAST/LAT based on the projected sensitivity 
for a one year observation. 

Key words. Acceleration of particles, ISM: cosmic rays, Galaxies: general, Galaxies: starburst, Galaxies: individual: M 82, 
Radiation mechanisms: non-thermal, Stars: formation 



1. Introduction 

In the central regions of starburst galaxies (SBGs) high level 
of star formation (SF) activity powers the emission of radiation 
directly by supernovae (SN) and massive stars, and indirectly 
by SN shock heating of interstellar gas and dust, as well as from 
radiative processes involving electrons and protons that are ac- 
celerated by SN shocks. Kinetic energy of the non-thermally 
distributed electrons and protons is partly channeled into radio 
to the high-energy 7-ray emission. In particular, radiation in 
the very-high-energy (VHE: >100 GeV) 7-ray region is emit- 
ted predominantly in 7r° decay, following the production of the 
pions in energetic proton interactions with protons in the gas, 
and by relativistic electrons traversing (relatively) strong mag- 
netic field and scattering off the intense radiation field in the 
galactic disk. 

The level of high-energy phenomena in intensely active 
galaxies makes AGN obvious targets of searches for VHE 
emission, and indeed some 15 blazars have already been de- 
tected by H.E.S.S. and MAGIC (e.g., Persic & De Angelis 
2008). Various issues in the study of non-thermal phenomena in 
'normal' and starburst galaxies are no less interesting and im- 
portant. Valuable insight can be gained on the origin and prop- 
agation modes of relativistic electrons and protons by prob- 
ing their radiation yields across the electromagnetic spectrum. 
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Knowledge of the low energy densities of these particles is im- 
portant also for quantifying their impact on interstellar gas. 

At a distance d=3.6 Mpc (Freedman et al. 1994), the 
SB galaxy M 82 is a good candidate for an initial search for 
TeV emission from a galaxy whose emission is not AGN- 
dominated. Its far-infrared (FIR) and X-ray luminosities imply 
high SF rates (e.g., Persic & Rephaeli 2007), and a correspond- 
ingly high SN rate of 0.3 yr^ 1 (Rieke et al. 1980). In most 
respects, levels of activity in M82 surpass those in the other 
nearby starburst galaxy, NGC 253, making it a very good tar- 
get for TeV observations, also in view of its large reservoir of 
diffuse gas (~lO 9 Af ; Casasola et al. 2004). 

In order to realistically estimate non-thermal emission from 
cosmic -ray electrons and protons, we start with the initial spec- 
trum that characterizes their momentum distribution at the ac- 
celeration sites, and follow the evolution of the spectrum as 
the particles diffuse and are convected from the inner SB re- 
gion to the outer disk (and halo). All the relevant leptonic and 
hadronic interactions are included in a numerical program we 
have developed for this purpose. Using this program we can de- 
termine the particle spatial and energy distributions, and their 
respective radiative yield across the galaxy. For a more intu- 
itive approach, and in order to compare with previous semi- 
quantitative estimates, we also use a simple analytical model 
of the diffuse hadronic emission of M 82, that ties the predicted 
7-ray flux to the measured on-site synchrotron radio flux by 
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assuming energetic particles and magnetic field to be in energy 
equipartition. 

In Section 2 we briefly review the basic features of pre- 
dicted particle source spectra, the relation between the cosmic- 
ray electron and proton components, and key aspects of their 
energy-loss processes and propagation in galaxies. We then 
summarize the main features of our numerical program which 
is used to determine the particle spectro-spatial distribution. 
In Section 3 we present a specific estimate of the hadronic 
emission from M 82 based on our approximate analytical treat- 
ment. Results of our detailed numerical analysis are presented 
in Section 4, and discussed in Section 5. 



2. Particle spectra and energy-loss processes 

2. 1. Particle source spectra 

SN shock wave sweeps and compresses circumstellar material 
with a compression ratio s=(7 c +l)/[7c _ 1+2 (f3/M) 2 ], where 
j c =Cp /Cv=5/3 is the specific heat ratio, (3=(c s /Va) is the 
speed of sound (c s ) in units of the Alfven velocity (Va), and 
M = ui/Va is the Alfvenic Mach number of the shock, with 
Ui denoting the upstream velocity of the gas. Particle accel- 
eration by the shock modifies the momentum distribution to a 
power-law whose index depends on Al. When particles with 
initial momentum p are overtaken by the shock, the result- 
ing downstream phase space distribution function is f(p) oc 
p~" (for p> Po ) where ?7=3s/(s-l)=3A/ 2 (7+l)/2(M 2 -l) 
(Krymsky 1977, Axford et al. 1977, Bell 1978, and Blandford 
& Ostriker 1978). The ensuing number density of accelerated 
particles per momentum interval dp is N(p) oc f(p)p 2 oc p~ q 
with q=r\—2. It can be readily seen that in the relativistic 
regime the particle density per unit energy interval has the same 
index, q. The strong shock limit has r)=4, i.e. q=2, a well- 
known result for the first-order Fermi acceleration mechanism 
(e.g., Protheroe & Clay 2004). 

Particles can be accelerated to very high energies (>10 14 
eV) even by the non-relativistic shocks in SN remnants 
(Aharonian et al. 2007), which are the accelerators of relevance 
to us here. The steepening of the particle power-law distribu- 
tion to q > 2 (due to energy losses outside the source region) 
implies that - for our purposes here - the specific value of the 
high-energy cutoff is irrelevant. For non-relativistic shock ac- 
celeration the lowest kinetic energy can be assumed to be the 
thermal energy of the gas, which is the source of accelerated 
particles. 

It is perhaps useful to outline the following simple - 
if physically uncertain - prediction (Schlickeiser 2002; Bell 
1978) for the proton to electron (p/e) ratio: Assume that elec- 
trons and protons are accelerated out of a thermal (Trj~/ew 
keV), electrically neutral plasma reservoir to high (T>Tq) ki- 
netic energies according to the same differential power-law 
momentum spectrum, Ni(p)—Nifip~ q with i=e,p, maintaining 
charge neutrality. 



1 It may be argued that charge neutrality does not necessarily re- 
quire equal numbers of positive and negative charges to be accelerated. 
Very minute current flows in the thermal plasma can easily compen- 



Under these assumptions, the p/e number and energy den- 
sity ratios in the source region can be analytically calculated. 
As function of kinetic energy T—(j~ l)mc 2 , where 7 is the 
particle Lorentz factor, the p/e number density ratio is found to 
be 



N P (T) 
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with m c and m p the electron and proton masses. This ratio 
assumes the following limiting forms: 
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In particular, the ratio reaches its maximum value, 
(m p /m e Y q ~ 1 ^ 2 (for q>l), over most of the range of 
particle energies. 

Similarly, the p/e energy density ratio, 
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can be expressed in terms of the spectral index q, 

n(q) = 
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For example, if q=2.3, then U p /U c ~15, and N p /N c ~l3 x 10 2 
at T> > 1 Ge V. These are approximately the values measured 
for Galactic cosmic rays (e.g., Schlickeiser 2002). 

The spectrum of relativistic electrons is obviously more 
easily measured due to their much more efficient radiative 
losses. Most readily observed is electron synchrotron radio 
emission, which we will use to infer N c (-f) in the region where 
the radio emission is observed. To relate this quantity to the 
source spectrum we need to solve the kinetic equation describ- 
ing the propagation mode and energy losses incurred by elec- 
trons as they move out from their SN shock regions. Using the 
above predicted proton to electron ratio at their common re- 
gion, and accounting for proton propagation mode and energy 
losses, we can then deduce the proton spectrum from the elec- 
tron spectrum. Clearly, even if the propagation modes of pro- 
tons to electrons are similar, their losses are very different; thus, 
N p should be inferred from particle number conservation. 

2.2. Steady-state particle spectra 

The enhanced SF activity in a starburst galaxy lasts for O(10 8 ) 
yr, which is also the timescale over which global SN shock 



sate for any charge imbalance in the accelerated particles. Indeed, be- 
cause electrons have much higher losses and different diffusion prop- 
erties, some degree of charge neutrality could develop during propa- 
gation away from the acceleration sites. 
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acceleration is prevalent. Estimated times for propagation out 
of the disk and inner galactic halo are about an order of mag- 
nitude lower for low energy particles. High energy particles 
lose energy on much shorter timescales. Thus, particle spectra 
evolve considerably as they propagate out of their acceleration 
sites. As long as typical acceleration times are comparable to or 
shorter than loss times for sufficiently high energies, a steady 
state can be attained, with particle spectra which are obviously 
steeper than their respective source spectra. Strictly speaking, 
we assume that steady state is attained and proceed to solve 
the kinetic equation for iVj(7, R, z), where R and z are the 2D 
spatial radius and the coordinate perpendicular to the galactic 
plane, respectively. 

A full calculation of the particle steady state spectra is nec- 
essarily extensive as it must include all the important energy 
loss mechanisms and modes of propagation. We have employed 
the numerical code of Arieli & Rephaeli (Arieli & Rephaeli, 
in preparation), which is based on a modified version of the 
GALPROP code (Moskalenko & Strong 1998, Moskalenko et 
al. 2003), that has been developed for this purpose. The code 
solves the exact Fokker-Planck transport diffusion-convection 
equation (e.g., Lerche & Schlickeiser 1982) in 3D with given 
source distributions and boundary conditions for electrons and 
protons. Evolution of the particle energy and spatial distribu- 
tion function is based on diffusion, and convection by galactic 
wind. (We note that some of the additional processes included 
in the original code, such as diffusive re-acceleration, nuclear 
fragmentation, and radioactive decay, are not relevant for our 
purposes here.) Logarithmic coordinates are used for high spa- 
tial resolution. 

The main focus of this first phase of our work is the estima- 
tion of the high energy flux of M 82. An approximate estimate 
of the high energy 7-ray emission can be obtained by a simple 
and intuitive calculation that is based on sampling the electron 
energy spectrum from the main observable - the radio syn- 
chrotron flux. However, because the electron spectrum needs 
to be extrapolated to energies much higher than the ~ 1 — 10 
GeV range directly inferred from radio measurements, realistic 
estimates of the TeV emission necessitates a detailed (numeri- 
cal) treatment. To illustrate this, and to provide insight on the 
results from the more accurate numerical calculation, we first 
outline the intuitive calculation. 



2.3.1. Gauging particle spectra by radio emission 

The electron population consists mostly of directly acceler- 
ated (primary) electrons, which can be represented by a single 
power-law (in terms of the Lorentz factor) 



iV e ( 7 ) = N e , 0l ~ q 



(5) 



for 7i<7<72- The electron synchrotron flux from a spherical 
region of radius r s with magnetic field B, located at a distance 
d, is 



3 

/„ = 5.7 x 1Q- 22 -^N c , a(q) B^r x 



4 x 10 6 



erg/ (s cm 2 Hz) 



(6) 



where a(q) is defined and tabulated in, e.g., Tucker (1975). 

Setting „=lGHz and ^ (^)- 3 (^) 2 (^f 
from Eq.© the normalization of the electron spectrum is 



iVe.o = 5.7 x 10~ 15 ^a^r 1 ^ — 250— cm" 



(7) 



A second relation is required to separately estimate N c and 
B. This is provided by assuming that the energy density is 
equipartitioned between particles (electrons and protons) and 
the magnetic field, 



B 2 



(8) 



The electron energy density is U e =N ei om e c 2 f^ 2 7 1 ~ 9 d 7 , 
where 7l ~100 corresponds to the energy below which 
Coulomb losses exceed synchrotron losses (the detailed num- 
ber depends on the ambient, mostly thermal, density and mag- 
netic field; e.g., Rephaeli 1979). For q>2 and 7 2>>7i, we 
get U c ~N e flm c c 2 ~f 2 ~ q / (q— 2). Substituting the expression for 
JV Ci o from Eq.(|7]i, we have 



U c = 3.0 x 10~ 22 2502 tp- 



7i 



? +2 



B 2 erg cm . 



(q - 2) a(q) 

Equation (O, together with Eqs.OJ and (8j, then leads to 



(9) 
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G. (10) 



2.3. Energy-loss processes 

At energies below few hundred MeV, electrons lose energy 
mostly by Coulomb interactions with gas particles, leading to 
ionization of neutral and charged ions, and electronic excita- 
tions in fully ionized gas. At higher energies the dominant en- 
ergy losses are via synchrotron emission, and Compton scatter- 
ing by the FIR and optical radiation fields. These well known 
processes need no elaboration; the level of the ensuing emis- 
sion depends on the mean strength of the magnetic field, B, 
and the energy density of the radiation fields, which are speci- 
fied below. 



From the above relations it can be readily seen that if the 
flux, source size, distance, and the spectral index of the emit- 
ting electron population (q) are known, then the source mag- 
netic field - and hence the particle energy density - can be 
evaluated. For example, if q=2.5 (which corresponds to a typi- 
cal galactic radio spectral index of 0.75, and entails a=0.0852 
and k~8) and i/j = 1, we obtain _B cq =61^G, ?7 C =10.25 eV 
cm" 3 , C7p=82 eV cm" 3 , A^ =10- 4 cm" 3 . 



2.3.2. FIR radiation and Compton X-ray emission 

In intensely star-forming regions the stellar IMF is top-heavy 
(e.g., Mayya et al. 2007), so massive stars are proportionally 
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more abundant there than in more typical Galactic environ- 
ments. During most of their lifetimes massive stars are em- 
bedded in the dense, highly extinguished regions (giant molec- 
ular clouds) where they were born (e.g., Silva et al. 1998). 
The radiation emitted by these stars peaks in the UV, which 
is very efficiently absorbed by the embedding dusty clouds and 
is re-emitted in the FIR. Therefore the local radiation field is 
most intense in the FIR region, where it can be described as an 
isotropic, diluted, modified blackbody with temperature Td and 
spatial dilution factor Cdii: 

1 e 2 e 

with 0<cr<2, for which we choose a—\ (e.g., Goldshmidt & 
Rephaeli 1995), and eo corresponding to v=2 x 10 12 Hz (Yun & 
Carilli 2002, and references therein). For a given Ta, the energy 
density of the warm dust, E/FiR(Cdii), can then be computed; 
comparison with the observed value, U§j^ — Lpj| l /(7rr 2 c), 
where r s is the radius of the accelerating source, yields Can. 

When electrons whose energy distribution is as in Eq.© 
interact with this FIR radiation, the spectral emissivity of the 
isotropically Compton-scattered radiation is (e.g., Rybicki & 
Lightman 1979) 

He) = C da F( q) ^!Y7 2) (kT d )^ e -V (12) 
7T Ti c 

where F( g ) S g ^^^ r(^) C(*£*). with T and C 
denoting the Gamma function and the Riemann zeta function. 
For GeV electrons and typical dust temperatures, most of the 
emission is in the X-ray region. 

2.3.3. 7-ray emission from n decay 

Proton energy losses at low energies are also dominated by 
Coulomb interactions with gas particles. At energies above 
pion masses (—140 MeV), the main energy-loss process is by 
interactions with ambient protons; the yield from this process 
are neutral (tt°) and charged (ir ) pions. Neutral pions decay 
into photons, while the decays of 7r yield relativistic e + +e~ 
and neutrinos. 

The integral spectral emissivity from ir° decay is: 

9 V {> ITeV) = ft,(=y 3 ~" (^Ir 1 (H atom)- 1 S ~\U) 
x lev cm- 3 

(Drury et al. 1994), where r}—q+2. This emissivity is very 
strongly dependent on the proton power-law index; e.g., for 
momentum distribution indices 77=4.1, 4.3, 4.5, 4.7, we deduce 
(using tabulated values from Drury et al. 1994) that the >100 
GeV emissivity, g{> 100 GeV), equals 1.3 X 10~ 16 , 4.2x 10~ 17 , 
9.5xl0~ 18 , 1.9xl0" 18 (ergcirr 3 )- 1 (H atom) -1 s~\ respec- 
tively. 

The integrated 7-ray photon luminosity of a source with gas 
density n and proton energy density U p in a volume V is then 

L(>e) = I g(>e)nU p dV s- 1 (14) 
Jv 

where both g(> e) and U p depend on the relativistic protons' 
spectral slope, q. 



3. VHE emission from M82: approximate 
treatment 

In this section we estimate the level of VHE emission from 
hadronic interactions in M 82, the most prominent nearby (d = 
3.6 Mpc) starburst galaxy. The radio emission of M82 is ob- 
served to be power-law from i/ m i„=22.25 MHz to ^ max =41 
GHz with a flux of 10 Jy at 1 GHz, and the spectrum steepen- 
ing from a ~ 0.71 (Carlstrom & Kronberg 1991) in the central 
starburst (defined here as a region with a radius of 300 pc and 
height of 200 pc; e.g., Volk et al. 1996, Mayya et al. 2006) to 
a ~ 1 (Klein et al. 1988, Seaquist & Odegard 1990) in the 
outer disk. The high-frequency radio data from the central disk 
region are shown in Figure 1 . While there is some curvature in 
the radio spectrum, because we are mostly interested here in the 
high-energy end of the electron spectrum, this curvature is ig- 
nored and we fit the spectrum by a single power-law model. We 
infer that the synchrotron radio emission is from a population 
of relativistic electrons whose energy spectrum - in the region 
sampled by the observed radio emission - is a power-law with 
index g~2.42. 

Using the implied value of q, from Eqs.© and (|4|i we 
obtain k~10 and jV„/iVJ,_, , ~ 2xl0 2 . Estimates of the 

VI HlGcV 

Galactic values of these energy densities are J7 P ~0.3 eV cm -3 
and £/ — 0.03 eV cm~ 3 , hence k— 10 (see Schlickeiser 2002), 
and N p /N c r^lQ 2 at 1 GeV. The values of the two ratios are then 
quite similar in the two galaxies, in spite of their very different 
environments. 

This similarity can be understood by the following argu- 
ment. Typical SNR ages (hence, the duration of the accelera- 
tion processes), tsnr^ 2x10 4 yr, are much shorter than the 
characteristic energy-loss timescales of both electrons and pro- 
tons in the remnant (even for the strongest magnetic fields in 
SNRs allowed by current IACT data, B = O(10 2 ) ^G; e.g., 
Berezhko & Volk 2006), r B ~l Gyr and 0.5 Myr for GeV pro- 
tons and electrons, respectively. Because of this, and the scal- 
ing of acceleration processes and propagation by rigidity, pro- 
tons and electrons in the central disk region are expected to 
retain approximate similarity of their spectra at low to moder- 
ate energies (which dominate the particle distributions). Thus, 
the spectra have roughly the same spectral index both at injec- 
tion ((fr~2) and, after injection, while diffusing in the central 
disk with a diffusion coefficient ocp b with fo— 0.5 (q—qi+b^2.5) 
(e.g., Biisching et al. 2001). This will hold both in the very ac- 
tively star-forming central region of M 82 and in the relatively 
quiescent Galactic environment. 

Using g=2.42 (that entails a~0.09 and n=U p /U c ~10), 
Eqs.dTOb. (|7]i, (0, and the above-quoted parameters ap- 
propriate for M82 (which entail ip=7.2), we derive: 
73-106 fiG, V Oj0 ~1.4xl0- 4 cm" 3 , C/ c ~25 eV cm" 3 , and 
f/ p ~250eVcm-' 3 . 

3. 1. X-ray emission 

Nonthermal X-ray emission is mostly by Compton scattering 
of the electrons by the FIR radiation field. To estimate the 
level of this emission, we first determine the energy density of 
the FIR field, whose dilution constant, Can, which appears in 



5 



Eq.dTTb. is evaluated as follows. Using Td=45°K (Carlstrom 
& Rronberg 1991) in Eq.dTTb. the warm dust energy density is 
C^fir=5-55 Cdii 10" 8 erg cm -3 . On the other hand, the mea- 
sured flux is /fir— 6.5 x 10~ 8 erg cm~ 2 s~\ so that the obser- 
vationally deduced value is C/^=L°^/(7rr s 2 c)~1.26x 10~ 9 . 
Equating U m to U$ s yields C di i=0.023. (This value is slightly 
reduced if one corrects the observed FIR luminosity of M 82 for 
cirrus emission.) 

From Eq.(fT2l the differential Compton flux from the central 
SB of M82 is / c (e)~4x lO" 13 ^ 71 erg cm" 2 s" 1 keV" 1 . 

The emission is clearly in the X-ray range; e.g., electrons 
with 7 as low as 150 and 350 boost a 60/^m photon (repre- 
sentative of the FIR field) to, respectively, 2 keV and 10 keV. 
Given their steep energy distribution, only electrons of rela- 
tively low energy will be effective in generating a substan- 
tial Compton flux. Thus, electrons with 100<7<1000 scatter- 
ing off the 45 °K blackbody photons will boost their energies 
mostly to the 0.4-80 keV band. In particular, the predicted 2-10 
keV Compton flux is / c (2 - 10 keV)~10~ 12 erg cm" 2 s _1 . 

There is yet no unequivocal evidence for nonthermal X-ray 
emission in M 82. Chandra observations show that some of the 
2-10 keV emission in M 82 is diffuse, emanating from a region 
which roughly overlaps with the central starburst (Griffiths et 
al. 2000). The 2-10 flux of this diffuse emission is 1.4 x 10~ 12 
erg s , some 5% of the total meassured 2-10 keV emission 
(e.g., Cappi et al. 1999; Rephaeli & Gruber 2002). Given that 
some of this diffuse component may be thermal, as suggested 
by the presence of a substantial 6.7 keV Fe-K emission, the ob- 
served diffuse emission is a strong upper limit on our predicted 
lower Compton contribution. 

3.2. VHE ~{-ray emission 

The value q=2A deduced for the central SB region corre- 
sponds to r;=4.4, for which we determine the appropriate value 
of g from Table 1 of Drury et al. (1994), and obtain g{> 
e)=8.1xlO- 19 (e/TeV)- 14 (erg cm- 3 )- 1 atom- 1 s" 1 , which 
holds down to GeV energies (e.g., Torres 2004; Domingo- 
Santamaria 2006). 

A rough estimate of the VHE emission from M 82 can be 
made if it assumed that this emission comes mostly from 7r° 
decay following interactions of energetic protons with ambi- 
ent protons in the SB region. In this approximate treatment we 
also ignore contribution to the VHE emission by electron radia- 
tive processes. Assuming that acceleration occurs mostly in the 
central SB (of radius r s =300 pc), we thus take the proton en- 
ergy density to be U P (R) ~ 250(i?/r s )~ 2 eV cm -3 for R>r s , 
and consider separately the emission from the central SB and 
the disk: 

(i) In the central SB region (R<r s ) the hydrogen mass 
is mostly molecular, Mh 2 — 2xl0 8 Af Q (Weiss et al. 
2001). The estimated VHE photon luminosity is then 
L(> 100GeV)~2xlO 39 s _1 , which corresponds to 
/(> 100 GeV)~1.3 x 10- 12 cm" 2 s" 1 . 

(ii) The disk region external to the central SB, R>r s , has 
a flat thin disk gas distribution, T,(R) = S(0)e _i?/flD , 
with S(0)~7.5xl0 22 cm -2 the effective central density 



and Rp^O.82 kpc, the radial lengthscale (Mayya et al. 
2006). The total gas mass, M gas ~2.5xlO 9 M0, results from 
M H i~0.7xl0 9 Ai© and M H2 ~1.8xlO 9 M (Casasola et al. 
2004). The energetic proton energy density is assumed to have 
the form U P (R)=200 (R/r s )~ 2 eV cm" 3 . From 

i-oo -x 

L(>e) ~ 2irR 2 D U p Y,{0)g{> e) / dx s" 1 , 

^0.366 x 

we then obtain L(> 100 GeV)~1.9xl0 40 s -1 , that corre- 
sponds to /(> 100 GeV)~1.2 x 10~ n cm" 2 s -1 . 

This approximate treatment then yields a total photon lu- 
minosity, ~2.1x 10 40 s _1 , which translates to a total flux /(> 
100 GeV) ~ 1.3 x 10~ n cm~ 2 s -1 . 

Before accepting the above estimate as physically meaning- 
ful, we need to ascertain that energetic protons do indeed leave 
their SN acceleration sites, and that VHE photons propagate 
freely through the galactic disk. Indeed, even in large, R ~ 10 
pc, dense molecular clouds with typical ambient proton den- 
sity of order 10 2 cm~ 3 , the optical depth to pp interactions is 
very low, T pp ~c7 pp n p /~3 x 10 -4 , where the cross section for 
pp interactions is <t pp ~10 2 mb=10 -25 cm 2 (e.g., Eidelman et 
al. 2004). This estimate is in obvious agreement with the fact 
that energetic protons at the relevant energy range are actually 
detected on Earth. 

Absorption of VHE photons is also negligible; the optical 
depth to pair production by interactions with ambient photons, 
T 11 —n 1 a 11 l (where n 1 is the target photon number density, 
<r 77 is the 77 cross section, and / is the distance traveled), 
is also very small. For TeV photons, the peak cross section, 
er 77 = 1.7 x 10 -25 cm 2 , is for interaction with ~2/xm pho- 
tons. The 2/jm luminosity, L2^m=2xl0 43 erg s" 1 (e.g., Silva 
et al. 1998), implies a 2/jm photon density of ^0.5 cm~ 3 (with 
a disk radius ?- gal ^7 kpc), i.e. t 77 ~10 -3 . We conclude that ab- 
sorption of VHE photons in - and indeed also along the l.o.s. 
to - M 82 is negligible. 

4. Particle and radiation spectra in M82: full 
numerical treatment 

As specified above, the central SB region (which will also re- 
ferred to as the source region) with a radius of 300 pc and 
height of 200 pc (e.g., Volk et al. 1996), is identified as the 
main site of particle acceleration. The superposed particle spec- 
tra from all SN shock regions yields the power-law distribu- 
tions with index q=2 (as mentioned in sect.2.1). The theo- 
retically predicted N p /N e ratio is more likely to be valid in 
this source region, as is also the assumption of equipartition. 
Accordingly, we infer N p in the source region from N , which 
itself is deduced from radio measurements in the central disk, 
as explained below. 

Particles are assumed to be convected out of their sources 
by the observed wind with velocity of 600 km s _1 (Strickland 
et al. 1997). The convection velocity is assumed to increase lin- 
early with distance from the disk plane. This assumption is con- 
sistent with cosmic -ray MHD wind models (e.g. Zirakashvili et 
al. 1996). The diffusion coefficient is taken as a power-law in 
rigidity D xx oc (p/po) s , where po is a scale rigidity with a typ- 
ical index of 0.5. Particle energy losses and propagation outside 




Fig. 1. Radio measurements (from Klein et al. 1988) and best- 
fit power-law spectrum to the data above 0.75 GHz. 



Fig. 2. Steady state primary electron (solid line) and proton 
(dashed line) spectral density distributions in the central disk 
region. 



the source region are followed in detail using the code of Arieli 
& Rephaeli. To do so we need the neutral and ionized gas den- 
sities and their profile, as well as the spatial variation of the 
mean strength of the magnetic field. 

We use the neutral, molecular, ionized gas, and FIR radia- 
tion field parameters specified above. Additionally, we assume 
magnetic flux freezing in the (highly conductive) ionized gas, 
so that the local magnetic field strength is related to the ionized 
gas density, nmi, using the scaling Bctn^ u (e.g., Rephaeli 
1988). Measurements of radio emission from the central disk 
(Klein et al. 1988), shown in Figure 1 together with our best fit 
power-law with index of 0.71 ± 0.02 (which is consistent with 
the results of Klein et al. 1988, and Carlstrom et al. 1993), pro- 
vide the first relation between the electron density in this region 
and the magnetic field. Note that we ignore here the flattening 
indicated in the radio spectrum at ^<300 MHz (which is pos- 
sibly due to the increased relative strength of Coulomb losses 
at low electron energies) since we are primarily interested in 
the high energy behavior of the electron spectrum. The second 
required condition to fully determine particle spectra and the 
magnetic field is the assumption of equipartition in the source 
region. Due to the implicit dependences in the expression for 
the synchrotron flux, this condition is implemented iteratively 
to solve for N c , N p , and Bq. With the above parameters the 
central value of the magnetic field is Bq ~ 180 fiG. 

A radio index of ^0.7 implies that the electron spectrum in 
the central disk steepens appreciably to q c ^2A (from q = 2 in 
the source region). The steady-state electron and proton spec- 
tra in the SB region are shown in Fig. 2. At low energies both 
spectra are flat; stronger electron losses at E»\ GeV result 
in a steeper spectrum than that of protons, with the electron 
spectrum characterized by q c< \^2.6 at very high energies. 

Electron emissions by bremsstrahlung and Compton scat- 
tering are shown in Fig. 3; also shown is 7-ray emission from 
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Fig. 3. Integral emission spectrum from the central disk region 
of M 82. Radiative yields are from electron Compton scatter- 
ing off the FIR radiation field (dotted line), electron bremm- 
strahlung off ambient protons (dashed line), ir° decay follow- 
ing pp collisions (dashed-dotted line), and their sum (solid 
line). 

7T° decay (following pp collisions). As expected, the losses due 
to bremsstrahlung dominate the lower energy regime, whereas 
losses due to ir° decay dominate at higher energies. (While syn- 
chrotron emission extends to the X-ray region, it is negligible 
at much higher energies.) Our main interest here is the VHE 
emission, which - as anticipated - is mainly from the latter 
process. Specifically, we calculate the total flux from M82 to 
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be /(> 100GcV)~ 2.1 xl(T 12 cmr 2 s" 1 . The flux at e>50 
GeV is a factor of ^3 higher. 

While we are mostly interested here in the VHE photon 
emission, we note that the related neutrino flux (7r ± eventually 
produce e ± + v e + v c + v„ + ) at energies higher than 100 
GeV is -0.3/(>100GcV). 

The predicted VHE 7-ray flux is a factor ~6 lower than 
that estimated in the approximate treatment in sect. 3. 2. Given 
the nature of the approximations made, and the fact that spatial 
dependences of the particle densities and magnetic field were 
altogether ignored, this level of discrepancy is not surprising. 
The two treatments differ most in the description of the particle 
spatial profiles outside the source region. In the approximate 
treatment the impact of energy losses and propagation mode 
of the protons are not explicitly accounted for; this results in 
the unrealistically high relative contribution of the main disk 
(exterior to the SB region) to the TeV emission. Clearly, even a 
small degree of steepening of the proton spectrum in this region 
results in a significantly lower TeV emissivity (see, e.g., Table 
1 of Drury et al. 1994) and hence flux. 

5. Discussion 

Our main objective has been an improved estimate of the VHE 
emission from the nearby northern-sky starburst galaxy M 82. 
This has been obtained first by an approximate calculation 
of the hadronic emission alone, and then a more comprehen- 
sive and detailed treatment which follows the evolution of the 
steady state electron and proton spectra and their radiative 
yields, taking into account pertinent energy losses and prop- 
agation modes. (Note though that because VHE emission is 
largely produced in the central disk, details of our modeling 
of the diffusion coefficient and convection velocity are not that 
important in determining the particle spatial profile over this 
relatively small region.) Of course, only the full, more physi- 
cally meaningful and realistic treatment provides a consistent 
basis for predicting the full particle and radiation spectra. The 
fact that we account for all known energy loss processes sub- 
stantiates our expectation that the electron spectrum can be ex- 
trapolated to energies much higher than the range directly sam- 
pled by radio measurements. The approximate calculation is 
included here in order to give intuitive insight on the main con- 
tribution to the VHE emission, and in order to relate our work 
to previous similar treatments. 

As we have already noted, the factor of ^6 discrepancy be- 
tween our two estimates of the total VHE flux is given the ap- 
proximate nature of the analytic calculation. This is so mainly 
because the bulk of the VHE emission comes from a relatively 
small central region of M 82 of ^1 kpc radius, essentially the 
SB, due to steep decline of the emission outside this source re- 
gion, which is described in the analytical treatment by x~ 1 e~ x . 
Over this region evolution of the particle spectra is moderate. 
Moreover, most of the VHE emission is from ir° decay, hence a 
detailed calculation of the electron spectrum is not required. Of 
course, the approximate calculation is totally inadequate in pre- 
dicting the lower energy emission, especially so over the larger 
disk region. 



Compared to previous estimates of 7-ray emission from 
M82 (e.g., Akyiiz et al. 1991; Pohl 1994; Volk et al. 1996), 
our detailed treatment is more realistic by virtue of being more 
tuned to a range of observables. These include the neutral and 
ionized gas densities and their estimated profiles across the cen- 
tral galactic disk, as well as the anticipated spatial variation 
of the magnetic field across this region. Most importantly, the 
basic normalization of the electron spectrum is directly deter- 
mined from the radio spectrum measured in the central disk 
region. We assume particle-field energy equipartition and an 
initial, theoretically predicted N p /N e ratio in the central re- 
gion, where these are more likely to be valid, rather than across 
the full disk as in previous work. Finally, we include all the 
dominant radiative electron and proton processes. The calcula- 
tion of Torres (2004), who estimated emission from the super- 
SB galaxy Arp220, and those of Paglione et al. (1996) and 
Domingo-Santamaria & Torres (2005), who estimated VHE 
emission from the southern-sky starburst galaxy NGC 253, al- 
though similar are still appreciably different from our numer- 
ical treatment. Given the substantial differences between our 
approach and those adopted in most previous treatments, it is 
not very meaningful to carry out a more detailed comparison of 
specific differences - e.g., in parameter values and spatial pro- 
files, including those of the electron and proton densities, and 
mean magnetic field in the central disk region - and we avoid 
doing so here. 

An assessment is needed of the reliability of our basic re- 
sult - the predicted level of VHE emission. Clearly, this flux 
depends linearly on 7V P , which in turn is related to N c . The 
electron density was deduced from the measured radio emis- 
sion in the central disk region. Therefore, the evolution of the 
electron spectrum from the latter region to the central disk re- 
gion is mostly determined by synchrotron losses in the high 
magnetic field. From Eqs. ^ and (HJ it follows that N p ocB 2 q 
(which is obvious in the limit U P »U C ). The uncertainty in the 
estimated level of TeV emission stems mostly from this steep 
dependence. It is unlikely that the field is appreciably higher 
than our estimate (Bq ~ 180 /iG), but it could possibly be 
lower. For example, a value lower by a factor of f, would re- 
duce N p and tt° decays by a factor C, 2 . However, for a given 
measured radio flux N e would have to be higher by the fac- 
tor £( 1+l ?)/ 2 , which for a mean electron index of ~ 2.4 (over 
the central disk) would result in a boost of the bremsstrahlung 
and Compton yields by nearly the same factor (~ £ 2 ). As can 
be seen in Figure 3, at the ~ 100 GeV region, the yield from 
the latter processes is only a factor ~ 2 lower than emission 
from 7T° decay, implying, e.g., that if the field were reduced to 
Bo— 100 fiG, the total VHE emission would not change signifi- 
cantly. Thus, the net effect of a lower field strength on the VHE 
emission is small. Of course, the ir° decay yield is linearly de- 
pendent also on the ambient proton density, which we inferred 
directly from the measured mass in the central disk. We con- 
clude that our predicted VHE flux is unlikely to be appreciably 
over-estimated. 



2 By a numerical coincidence in spite of the different approach, our 
(numerical) VHE flux estimate is very similar to that by Volk et al. 
1996, which was assumed to be emitted from the central starburst 
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Based on the above estimate, the predicted VHE flux of 
M 82 falls somewhat below the detection limit of current imag- 
ing air Cherenkov telescopes (IACTs). For example, MAGIC 
and VERITAS - which are located in the Northern hemisphere 
and can therefore observe M 82 - have a 5a sensitivity of, re- 
spectively, ^2xl0~ n and 10~ n cm~ 2 s _1 for the detection 
of emission above 100 GeV in 50 hours of observation (e.g., 
Bastieri et al. 2005). However, this sensitivity limit refers to 
a Crab-like spectrum, i.e. a power-law source whose differen- 
tial VHE (photon) flux has an index of 2.6, somewhat steeper 
than our deduced value of ^2.3. Accounting for this and for 
the inherent uncertainty (by a factor of ^ 2-3) in our lower 
flux estimate, we conclude that measurement of VHE emission 
from M 82 - or obtaining an interesting upper limit on it - may 
be feasible with MAGIC or VERITAS, if observed for more 
than ^500 hours. Prospects for detection are real with the up- 
coming MAGIC II telescope, whose sensitivity at >100 GeV 
is expected to be a factor of ^3 better than that of MAGIC I. 

Finally, we comment on the possibility of detecting M 82 
at >100MeV 7-rays with GLAST. Our numerical code yields 
an integrated flux /(>100 MeV)~l.l x 10" 8 an^s" 1 . This 
value matches the 5<r sensitivity of GLAST& Large Area 
Telescope (LAT) for a 1-yr scanning-mode operation (e.g., 
Dermer 2007). Thus, GLAST/LAT could possibly detect M 82 
during its first year of operation. 

Its proximity and strong star-formation activity make 
M82 the most promising non-AGN galaxy target for high- 
energy 7-ray detection with the current generation of detectors 
(GLAST/LAT, MAGIC, and VERITAS). The LMC is the only 
non-AGN galaxy detected at energies >100MeV (Sreekumar 
et al. 1992); weak upper limits have been set on the VHE 
emission from normal and SB galaxies, including NGC253 
(Aharonian et al. 2005c) and Arp220 (e.g., Albert et al. 2007). 
Detection of VHE 7-ray emission from nearby star-forming 
galaxies is clearly of much interest in view of the feasibility 
of testing cosmic -ray acceleration and propagation models. 

Acknowledgements. We thank the referee for useful re- 
marks on some aspects of particle acceleration by non- 
relativistic shocks. 
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